﻿namespace DiffLCS

open System
open System.Threading.Tasks
open System.Collections.Generic
open TaskHelper

// TODO: Make the classes static
// TODO: Refactor common code from the parallel algorithms

[<Flags>]
type Settings =
    | Sequential = 1
    | Parallel = 2
    | Auto = 3

type LCS() =
    static let processorCount = Environment.ProcessorCount

    // TODO: Find an optimal level, maybe expose these
    [<Literal>]
    static let qLengthAuto = 500
    [<Literal>]
    static let minIterationsAuto = 1000000

    static member private runParallel n m settings =
        match settings with
        | Settings.Auto -> processorCount > 1 && n * m >= minIterationsAuto
        | Settings.Parallel -> true
        | Settings.Sequential -> false
        | _ -> raise (new ArgumentException("Invalid settings argument"))

    /// Returns a LCS backtrack matrix for the strings "a" and "b". This is done sequential.
    static member internal MatrixSeq (a : string, b : string) =
        let (n, m) = (a.Length, b.Length)
        let C = Array2D.zeroCreate (n + 1) (m + 1)
    
        for i = 1 to n do
            for j = 1 to m do
                C.[i, j] <- if a.[i - 1] = b.[j - 1] then
                                C.[i - 1, j - 1] + 1
                            else
                                Math.Max(C.[i, j - 1], C.[i - 1, j])
        C

    /// Returns a LCS backtrack matrix for the lists "a" and "b". This is done sequential.
    static member internal MatrixSeq (a : IList<'a :> IEquatable<'a>>, b : IList<'a>)  =
        let (n, m) = (a.Count, b.Count)
        let C = Array2D.zeroCreate (n + 1) (m + 1)
    
        for i = 1 to n do
            for j = 1 to m do
                C.[i, j] <- if a.[i - 1].Equals(b.[j - 1]) then // MUCH faster than =
                                C.[i - 1, j - 1] + 1
                            else
                                Math.Max(C.[i, j - 1], C.[i - 1, j])
        C

    /// Returns a LCS backtrack matrix for the arrays "a" and "b". This is done sequential.
    static member internal MatrixSeq (a : ('a :> IEquatable<'a>)[], b : 'a[]) =
        let (n, m) = (a.Length, b.Length)
        let C = Array2D.zeroCreate (n + 1) (m + 1)
    
        for i = 1 to n do
            for j = 1 to m do
                C.[i, j] <- if a.[i - 1].Equals(b.[j - 1]) then
                                C.[i - 1, j - 1] + 1
                            else
                                Math.Max(C.[i, j - 1], C.[i - 1, j])
        C

    /// Returns a LCS backtrack matrix for the strings "a" and "b". This is done in parallel.
    static member internal MatrixPar (a : string, b : string, qLength) =
        let (n, m) = (a.Length, b.Length)

        let C = Array2D.zeroCreate (n + 1) (m + 1)

        if a = "" || b = "" then
            C
        else
            // Count of how many parts (qx * qy) the computation is divided into
            let qx =  Math.Max(n / qLength, 1)
            let qy = Math.Max(m / qLength, 1)

            // Returns a 4-tuple representing the part q_x,y
            // TODO: Optimize so length is n / qx and height is m / qy.
            // This means qLength will be a approximation.
            let qCoord x y =
                let i = x * qLength + 1
                let j = y * qLength + 1
                let w = if x = qx - 1 then
                            n - x * qLength
                        else
                            Math.Min(qLength, n)
                let h = if y = qy - 1 then
                            m - y * qLength
                        else
                            Math.Min(qLength, m)

                (i, j, w, h)
    
            // Performs the actual computation for a specific part q_x,y
            let calculateLCSMatrix x y =
                let (i, j, w, h) = qCoord x y
                let (i1, j1) = (i + w - 1, j + h - 1)

                for k = i to i1 do
                    for l = j to j1 do
                        C.[k, l] <- if a.[k - 1] = b.[l - 1] then
                                        C.[k - 1, l - 1] + 1
                                    else
                                        Math.Max(C.[k, l - 1], C.[k - 1, l])

                ()

            let calculateFun x y =
                (fun () -> calculateLCSMatrix x y)

            // An array of tasks, each representing a computation of a part q_x,y
            let Q = Array2D.create qx qy (null :> Task)

            // The following few lines setups the computation of different parts,
            // see the definition of --> and >-> in TaskHelper

            Q.[0, 0] <- (task (calculateFun 0 0))

            for x = 1 to qx - 1 do
                Q.[x, 0] <- Q.[x - 1, 0] --> calculateFun x 0

            for y = 1 to qy - 1 do
                Q.[0, y] <- Q.[0, y - 1] --> calculateFun 0 y

            for x = 1 to qx - 1 do
                for y = 1 to qy - 1 do
                    Q.[x, y] <- [| Q.[x - 1, y]; Q.[x, y - 1] |] >-> calculateFun x y

            // Starts the first computation which the rest depends on
            Q.[0, 0].Start();
            // Wait for the last computation to finish
            Q.[qx - 1, qy - 1].Wait();

            C

    /// Returns a LCS backtrack matrix for the strings "a" and "b". This is done in parallel.
    static member internal LengthPar (a : string, b : string, qLength) =
        let (n, m) = (a.Length, b.Length)

        if a = "" || b = "" then
            0
        else
            // Count of how many parts (qx * qy) the computation is divided into
            let qx =  Math.Max(n / qLength, 1)
            let qy = Math.Max(m / qLength, 1)

            // Returns a 4-tuple representing the part q_x,y
            // TODO: Optimize so length is n / qx and height is m / qy.
            // This means qLength will be a approximation.
            let qCoord x y =
                let i = x * qLength + 1
                let j = y * qLength + 1
                let w = if x = qx - 1 then
                            n - x * qLength
                        else
                            Math.Min(qLength, n)
                let h = if y = qy - 1 then
                            m - y * qLength
                        else
                            Math.Min(qLength, m)

                (i, j, w, h)
    
            // Performs the actual computation for a specific part q_i,j
            let calculatePart x y (c00 : int[,]) (c10 : int[,]) (c01 : int[,]) =
                let (i, j, w, h) = qCoord x y

                let (w0, h0) = (Array2D.length1 c00, Array2D.length2 c00)
                printf "Mem: %d, " (GC.GetTotalMemory(true))
                let c11 = Array2D.zeroCreate w h

                c11.[0, 0] <- if a.[i - 1] = b.[j - 1] then
                                  c00.[w0 - 1, h0 - 1] + 1
                              else
                                  Math.Max(c10.[w0 - 1, 0], c01.[0, h0 - 1])

                for k = 1 to w - 1 do
                    c11.[k, 0] <- if a.[i + k - 1] = b.[j - 1] then
                                      c01.[k - 1, h0 - 1] + 1
                                  else
                                      Math.Max(c01.[k, h0 - 1], c11.[k - 1, 0])

                for l = 1 to h - 1 do
                    c11.[0, l] <- if a.[i - 1] = b.[j + l - 1] then
                                      c10.[w0 - 1, l - 1] + 1
                                  else
                                      Math.Max(c10.[w0 - 1, l], c11.[0, l - 1])

                for k = 1 to w - 1 do
                    for l = 1 to h - 1 do
                        c11.[k, l] <- if a.[i + k - 1] = b.[j + l - 1] then
                                          c11.[k - 1, l - 1] + 1
                                      else
                                          Math.Max(c11.[k, l - 1], c11.[k - 1, l])

                c11

            // The following few lines setups the computation of different parts,
            // see the definition of -->> and >->> in TaskHelper

            let rightZeroArr y =
                let (_, _, _, h) = qCoord 0 y

                Array2D.zeroCreate 1 h

            let bottomZeroArr x =
                let (_, _, w, _) = qCoord x 0

                Array2D.zeroCreate w 1

            let construct () =
                // An array of tasks, each representing a computation of a part q_x,y
                let Q = Array2D.create qx qy (null :> Task<int[,]>)

                Q.[0, 0] <- (taskf (fun () -> (calculatePart 0 0 (Array2D.zeroCreate 1 1) (rightZeroArr 0) (bottomZeroArr 0))))

                for i = 1 to qx - 1 do
                    Q.[i, 0] <- Q.[i - 1, 0] -->> (fun q -> calculatePart i 0 (bottomZeroArr (i - 1)) q (bottomZeroArr i))

                for j = 1 to qy - 1 do
                    Q.[0, j] <- Q.[0, j - 1] -->> (fun q -> calculatePart 0 j (rightZeroArr (j - 1)) (rightZeroArr j) q)

                for i = 1 to qx - 1 do
                    for j = 1 to qy - 1 do
                        Q.[i, j] <- [| Q.[i - 1, j]; Q.[i, j - 1] |] >->> (fun (qs : int[,][]) -> calculatePart i j Q.[i - 1, j - 1].Result qs.[0] qs.[1])

                (Q.[0, 0], Q.[qx - 1, qy - 1])

            let (first, last) = construct ()

            // Starts the first computation which the rest depends on
            first.Start();

            // Wait for the last computation to finish
            last.Wait();
            
            last.Result.[last.Result.GetLength(0) - 1, last.Result.GetLength(1) - 1]

    /// Returns the length of a longest common subsequence.
    /// Space: O(min{|a|, |b|})
    static member Length (a : string, b : string) =
        // If a is longer than b, then optimize memory by swapping.
        let (a, b) = if b.Length > a.Length then
                        (b, a)
                     else
                        (a, b)
        
        let (n, m) = (a.Length, b.Length)

        let C = Array2D.zeroCreate 2 (m + 1)

        for i = 1 to n do
            for j = 1 to m do
                C.[1, j] <- if a.[i - 1] = b.[j - 1] then
                                C.[0, j - 1] + 1
                            else
                                Math.Max(C.[1, j - 1], C.[0, j])

            for j = 1 to m do
                C.[0, j] <- C.[1, j]

        C.[1, m]

    static member internal Sequence (a : string, b : string, backtrack : int[,]) =
        let C = backtrack

        let rec visit i j =
            if i = 0 || j = 0 then
                []
            elif a.[i - 1] = b.[j - 1] then
                a.[i - 1] :: visit (i - 1) (j - 1)
            elif C.[i, j - 1] > C.[i - 1, j] then
                visit i (j - 1)
            else
                visit (i - 1) j

        visit a.Length b.Length |> List.rev
                                |> List.toString

    static member Sequence (a : string, b : string, settings : Settings) =
        let C = if LCS.runParallel a.Length b.Length settings then
                    LCS.MatrixPar(a, b, qLengthAuto)
                else
                    LCS.MatrixSeq(a, b)

        LCS.Sequence(a, b, C)
        

    static member Sequence (a : string, b : string) =
        LCS.Sequence(a, b, Settings.Auto)